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1.  INTRODUCTION 


The  accurate  simulation  of  the  flowfield  of  a  helicopter  rotor  is  still  one  of  the  most 
complex  and  challenging  problems  of  applied  aerodynamics.  This  is  true  in  spite  of  the 
availability  of  the  present  day  supercomputers  of  Cray-2  class  and  improved  numerical 
algorithms.  The  major  reason  for  this  is  that  the  flowfield  of  a  rotor  in  forward  flight  is 
highly  three-dimensional,  unsteady  and  viscous,  with  pockets  of  transonic  flow  near  the 
blade  tips  on  the  advancing  blades  and  regions  of  dynamic  stall  on  the  retreating  blades. 
In  addition,  the  blades  also  shed  complex  vortical  wakes;  the  concentrated  tip  vortices 
shed  by  these  blades  generally  pass  close  to  the  following  blades.  The  close  encounter 
of  these  (force  free)  concentrated  vortices  with  the  rotor  blades  induces  unsteady  load 
fluctuations  and  acoustic  noise.  Some  of  the  complexities  associated  with  such  a  flowfield 
are  delineated  in  a  schematic  of  a  helicopter  rotor  in  forward  flight  in  Fig.  1. 

The  spiralling  vortex  sheet  emanating  from  each  of  the  blades  of  the  rotor  has  a 
profound  influence  on  the  performance  of  the  helicopter.  It  not  only  alters  the  effective 
pitch  angle  of  the  blades  and  thus  the  airloads,  but  also  produces  highly  nonlinear 
interactions  of  the  vortex  with  the  rotor  flowfield.  It  is  possible  that  such  interactions 
might  produce  vortex  induced  boundary  layer  separation  which  results  in  a  sudden  loss 
of  lift  and  increase  in  drag.  An  accurate  simulation  of  the  rotor  flowfield,  therefore, 
must  consider  the  induced  effects  of  the  vortex  wake  including  the  effects  of  blade- vortex 
interactions  (BVI). 

Numerical  simulations  of  the  vortex  wakes  are  being  attempted  only  recently  after 
bigger  and  faster  computers  have  become  available  (Ref.  1).  These  investigations  have 
had  some  limited  success  to  date  in  simulating  the  formation  and  initial  roll  up  processes. 
On  the  other  hand,  much  of  the  progress  in  modeling  the  blade-vortex  interaction  has 
been  hampered  by  the  lack  of  development  of  theoretical  and/or  numerical  techniques 
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to  preserve  the  structure  of  the  concentrated  vortices,  when  convecting  in  the  flowfield, 
without  significant  diffusion.  The  study  of  blade- vortex  interaction  has  been  the  subject 
of  numerous  recent  research  papers  (Refs.  2-10).  These  studies  have  considered  different 
methods  of  preserving  the  concentrated  vortex  in  the  flowfield.  Among  these,  the  vortex¬ 
fitting  technique  (also  called  the  split  potential  formulation  in  potential  flow  methods) 
has  been  demonstrated  to  be  a  very  effective  method  in  preserving  the  vortex  even  when 
the  computational  grid  is  sparse.  This  method  has  also  proven  very  economical  (Ref.  5) 
compared  to  a  more  exact  formulation  used  in  Ref.  3. 

Accurate  simulations  of  the  flowfield  of  a  helicopter  rotor  in  hover  or  in  forward  flight 
require  an  accurate  description  of  the  vortex  structure  and  its  position  in  the  wake. 
The  Navier-Stokes  methods  are  the  only  methods  that  can  provide  this  information. 
To  compute  such  an  information  on  a  routine  basis  for  a  design  analysis  is  perhaps 
prohibitively  expensive  even  with  the  present  day  supercomputers.  This  is  the  reason 
that  simple  methods,  based  on  linear  aerodynamics,  have  been  developed  for  use  in 
industry  (see  for  example  Refs.  11-12). 

To  date,  many  advances  have  been  made  to  calculate  these  complex  flowfields  by  the 
use  of  simpler  equations  of  fluid  motion,  such  as  full  potential,  coupled  with  empirical 
wake  models  to  bring-in  the  influence  of  the  vortex  wake  (Ref.  8).  Such  methods  have 
been  specialized  for  each  blade  shape  and  planform  and  therefore  are  not  very  general. 
Nevertheless,  these  methods  have  been  useful  to  the  industry  and  are  being  used  in 
the  design  analysis  on  a  routine  basis.  There  is  no  doubt  that  the  weak  link  in  the 
above  methodology  is  the  wake  modeling  as  these  methods  cannot  handle  arbitrary 
configurations  and  high  supercritical  flow  conditions.  Alternative  to  the  prescribed  wake 
models  used  in  the  above  methods,  the  present  study  will  present  a  numerical  method  to 
capture  the  wake  and  its  influence  as  a  part  of  the  overall  flowfield  solution.  This  means 
that  the  present  method  does  not  use  any  ad  hoc  wake  models. 

It  should  be  pointed  out  that  the  concept  of  wake  capture  is  not  new  and  has,  in 
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fact,  been  utilized  by  several  researchers  to  investigate  the  hovering  rotor  problem.  In 
particular,  the  works  of  Refs.  13-17  are  to  be  mentioned  for  the  Euler  and  flow  potential 
methods.  These  studies  have  shown  promise  with  limited  success.  The  main  ingredient 
that  is  responsible  for  the  formation  of  tip  vortices,  namely  the  viscosity,  is  absent  in 
these  methods.  While  the  role  of  viscosity  is  clear,  its  importance  in  these  flows  is 
uncertain  at  this  stage.  This  will  be  investigated  in  the  present  study. 

The  purpose  of  this  study,  therefore,  is  to  devise  numerical  methods  for  the  solu¬ 
tion  of  the  Navier-Stokes  equations  to  accurately  calculate  the  flowfield  of  a  helicopter 
rotor,  including  the  effects  of  blade-vortex  interactions.  In  particular,  this  report  will 
attempt  to  demonstrate  the  ability  of  the  current  Navier-Stokes  methods  to  1)  calcu¬ 
late  the  unsteady,  three-dimensional,  blade-vortex  interaction  (BVI)  flowfield,  for  both 
parallel  and  oblique  blade- vortex  interactions  for  a  model  helicopter  rotor  problem  (see 
for  example  Refs.  18-19),  and  2)  calculate  the  flowfield  of  a  lifting  hovering  rotor  by 
capturing  the  wake  and  its  influence  as  a  part  of  overall  solution  method  corresponding 
to  the  experimental  conditions  of  a  model  rotor  test  (Ref.  20). 
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2.  GOVERNING  EQUATIONS  AND  NUMERICAL  SCHEME 


The  governing  partial  differential  equations  are  the  unsteady,  thin  layer  Navier- 
Stokes  equations.  In  this  study,  the  thin  layer  Navier-Stokes  code  recently  written  specif¬ 
ically  for  helicopter  rotor  applications  (Ref.  21)  will  provide  the  starting  point  for  further 
modifications  to  calculate  the  blade- vortex  interaction  flowfield  first  and  later  the  lifting 
hovering  rotor  flowfield.  In  order  to  assess  the  importance  of  viscous  effects  in  the  blade- 
vortex  interacting  flowfields,  both  the  Euler  and  the  thin  layer  Navier-Stokes  equations 
have  been  solved.  (The  Euler  calculations  are  performed  by  turning  off  the  viscous  terms 
in  the  Navier-Stokes  equations  set  and  appropriately  modifying  the  surface  boundary 
condition.)  The  mathematical  formulation  presented  in  the  following  paragraphs  is  for 
the  solution  of  thin  layer  Navier-stokes  equations  with  the  added  feature  of  introducing 
a  prescribed  vortex. 

For  generality,  the  equations  are  transformed  from  the  inertial  Cartesian  refer¬ 
ence  frame  ( x ,  y,  z,  t )  to  the  arbitrary  curvilinear  coordinate  frame  (£,  ij,  (,  r)  that 
moves  with  the  blade,  while  retaining  strong  conservation  law-form  to  capture  shock 
waves.  The  transformed  equations,  following  Ref.  22,  can  be  written  as 


dr(Q  -  Q0)  +  d((F  -  Fo)  +  6V(G  -  Go)  +  d((H  -  H0)  =  tRe~' d<S  (1) 


where 
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V(e+p)  -  T]tp 

W(e  +  p)  -  CtP 

and  e  =  0  for  Euler  equations,  and  e  =  1  for  thin-layer  Navier-Stokes  equations.  Here 

/v  /s 

Q  is  the  unknown  flowfieid  vector  and  Q0  is  the  solution  of  the  Elder  equations  for  a 
prescribed  line  vortex  aligned  with  the  uniform  free  stream  of  Mach  number,  Af^,  in  the 

/V  ^  XV  XV  XV  XV 

y-direction  (see  Fig.  2).  The  terms  F,  G,  H  and  F0,  G0 ,  H0  represent  the  appropriate 
flux  vectors  for  the  two  flows,  respectively.  The  hat,  represented  by  (  ),  in  Eq.  (1) 

denotes  quantities  scaled  by  the  Jacobian  J,  e.g.,  Q  =  J~*Q  etc.  The  contravariant 
velocity  components  U,  V  and  W  are  defined  as 

U  =  (t  +  (*u  4-  iyv  +  £zw 

V  =  T]t  +  rixu  +  T]yV  +  TJZW 

W  =  Ct  4-  (xu  4-  Cyv  4-  C,zw 

In  the  present  formulation,  £  lies  in  the  chordwise  or  wraparound  direction,  y  is  in  the 
spanwise  direction,  and  £  is  normal  to  the  blade  surface. 

XV 

The  viscous  flux  vector  5  is  used  in  the  limit  of  the  thin-layer  approximation 
and  is  given  by 


0 

Kiuc  +  K2(t 
KiV(  +  K2Cy 
KiW(  -I-  K2C,z 

l)-V)<+«72)/2)<) 
+  K2Kz 


(3) 
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where 

Ki  =  kci + <; + a 

#2  =  MC*UC  +  Cy*><  +  Cz^<)/3 

(4) 

-K3  =  C*u  +  Cyv  +  Cz™ 
y2  =  u2  +  v2  +  to 2 

In  this  equation,  a  is  the  speed  of  sound,  k  is  the  coefficient  of  thermal  conductivity  and  p 
is  the  viscosity  coefficient  computed  as  the  sum  of  laminar  and  turbulent  contributions. 
Sutherland’s  equation  is  used  for  laminar  viscosity  and  the  two-layer  algebraic  eddy 
viscosity  model  of  Baldwin  and  Lomax  (Ref.  23)  is  used  to  calculate  the  turbulent  eddy 
viscosity,  after  minor  modifications  to  account  for  the  moving  blade  surface.  Also,  J  is 
the  transformation  Jacobian,  whose  inverse  is  written  as 


J  1  =  +  X^Zr,  +  Xr,y<Zi  -  X^Zj,  -  X^Z^  -  X(T/,Zf 

The  velocity  components  u,v,w  and  the  pressure,  p,  are  related  to  the  total 
energy  per  unit  volume,  e,  through  the  equation  of  state  for  a  perfect  gas  by 


V  =  (7  “  l)(e  ~  ^(«2  +v2  +  w2))  (5) 

The  primitive  variables  of  Eq.  (1),  viz.,  the  density  p,  the  mass  fluxes  pu,pv,pw 
and  the  energy  per  unit  volume  e,  are  normalized  by  the  free  stream  reference  quantities 
and  the  pressure  p  by  7Poo-  The  reference  length  and  velocity  scales  are  the  chord  of  the 
blade,  C,  and  free  stream  speed  of  sound,  a <*,,  respectively.  The  details  of  the  metrics  of 
transformation  (&>£*> {TH,‘nz,T)y,Tlz)  and  (Ct, C*> Cy> Cz)  can  be  found  in  Ref.  22. 

In  the  above  equations  tx,v,  and  w  are  the  Cartesian  components  of  velocity  in 
the  inertial  coordinate  system  (x,y,  z,t).  In  the  present  formulation  Eq.  (1)  is  solved 
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in  the  inertial  frame  of  reference  with  the  boundary  conditions  applied  on  the  rotating 
blade.  The  inertial  coordinates  X  =  (a;,y,z,t)  are  related  to  the  blade-fixed  coordinates 
Xb  =  (x,y,z)  through  the  relation  given  by 


X{x,y,z)  =  R(t)Xh(x,y,z) 
t  =  f 

where  R(t )  is  the  rotational  matrix  (Ref.  24)  given  by 


(6) 


cos  fit  —  sin  f If  0 

R(t)  =  sin  fit  cos  fit  0  (7) 

0  0  1 


Here  fl  is  the  angular  velocity  of  the  rotor  and  fit  represents  the  azimuth  sweep  of  the 
rotor  blade. 


The  equations  set,  Eq.  (1),  is  solved  using  a  two-factor,  implicit,  finite- difference 
numerical  scheme  (Ref.  25),  that  uses  spatial  central-differencing  in  the  tj  and  £  directions 
and  upwind-differencing  in  the  £  direction  and  is  given  by 

[/  +  h8l(A+)n  +  -  hRe~16('J-iMnJ  -  Z>*|c] 

x  [/  +  h6'(A~)n  +  hSvBn  -  Dil,]  (A Qn  -  A $»)  = 

-  A t{6b[(F+)n  -  (F0+)n|  +  5'[(F-)n  -  (K)n)  (8) 

+  6v(Gn  -  Gn0)  4-  SdHn  -  H:)  -  Re-1 6^)} 

-  (De\v  +  De\()(Qn  -  &) 

where  h  =  At  for  first-order  time  accuracy.  Here  6  is  typically  a  three-point,  second- 
order-accurate,  central  difference  operator,  6  is  a  midpoint  operator  used  with  the  viscous 
terms,  and  the  operators  and  6 £  are  backward  and  forward  three-point  difference 
operators.  The  time  index  is  denoted  by  n  such  that  t  =  (nAf),  and  A Qn  =  Qn+1  —  Qn. 

A  “V 

The  flux  vector  F  has  been  split  into  F+  and  F~ ,  according  to  its  eigenvalues  (Ref.  26) 
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and  the  matrices  A±,  B ,  C ,  and  M  result  from  the  local  linearization  of  fluxes  about 
the  previous  time  level  (see  Ref.  22).  The  finite-difference  scheme  described  in  Eq.  (5) 
uses  flux-splitting  in  £  direction  and  central  differencing  in  the  77  and  £  directions.  As 
a  consequence,  numerical  disspation  terms  denoted  by  Di  and  De  are  employed  in  the 
77  and  £  directions,  and  axe  given  as  combinations  of  second-  and  fourth- differences  as 
follows  (Ref.  27). 


Del,  =  +  \V\y  +  \V\Z) 

D,|,  =  (A0J-1(kU  +  |7?|y  +  |7?U) 


£28 

£28 


\62p\ 

\(1+P)p\ 

\&p\ 

l(l  +  *2)pl 


8  +  e<84 


8  +  3e4£4 


(9) 


In  the  vicinity  of  shock  waves,  the  fourth-difference  terms  can  cause  oscillations,  so  it  is 
desirable  to  drop  these  terms  and  rely  only  on  the  second-difference  terms.  In  the  above 
expression,  p  is  the  nondimensionalized  fluid  pressure  and  e2  is  0(1)  while  e4  is  0(0.1). 
In  transonic  flow  the  term  (|?7|*  +  |?7 1 y  +  |?7|2)  is  an  estimate  to  the  spectral  radius  of 
the  matrix  B  for  variables  that  have  been  nondimensionalized  with  respect  to  the  sound 
speed,  since  (|q|z  +  ItIs  +  Mz)  does  not  contain  any  fluid  variable,  it  can  be  left  as  the 
coefficient  to  the  dissipation  operator  without  affecting  the  weak  conservation  form  of 
the  difference  equations.  The  implicit  second  difference  numerical  dissipation  operators 
were  chosen  to  keep  the  left  hand  side  factors  block  tridiagonal,  and  as  D,  works  on  8Q, 
accuracy  is  not  impaired. 


The  factored  operators  in  Eq.  (8)  are  solved  by  sweeping  in  the  £  direction  and 
inverting  tridiagonal  matrices  with  5x5  blocks  for  the  other  two  directions.  Currently,  a 
significant  part  of  the  computational  time  is  taken  to  form  the  plus  and  minus  Jacobian 
matrices  for  the  flux  vector  F  with  this  numerical  scheme.  However,  this  effort  has  been 
reduced  by  computing  and  A"  at  every  other  point  (in  both  77  and  (  directions) 
and  averaging  to  obtain  the  matrices  at  the  intermediate  points.  The  numerical  code  is 
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vectorized  for  the  Cray-2  supercomputer. 

In  writing  Eq.  (5)  it  is  assumed  that  .A*  ssA,  B0  ~  B ,  and  C0  rs  C  where 
the  Jacobian  matrices  ,  B0  and  C0  correspond  to  the  prescribed- vortex  flowfield.  In 
the  absence  of  vortex  interaction,  the  prescribed-disturbance  flowfield  reduces  to  a  free 
stream. 

In  using  this  numerical  algorithm  for  the  physical  problem  of  blade-vortex  in¬ 
teraction  presented  here,  a  note  of  caution  is  in  order.  As  noted  before,  the  numerical 
scheme  is  second-order  accurate  in  space.  Although  it  is  desirable  to  have  higher-order 
accuracy  in  space,  it  turns  out  that  the  present  algorithm  is  numerically  more  stable  if 
this  accuracy  is  decreased  slightly  from  second-order.  Also,  in  evaluating  the  right  hand 
side  of  Eq.  (5),  it  is  desirable  to  calculate  the  derivatives  of  individual  fluxes  first  and 
then  taking  a  difference.  This  is  particularly  true  for  the  fluxes  in  the  £  direction. 

The  above  numerical  algorithm  is  used  for  the  problem  of  blade-vortex  inter¬ 
action.  For  calculating  the  lifting  hovering  flowfields,  a  slightly  different  numerical  al¬ 
gorithm  is  used.  The  fundamental  difference  is  the  use  of  upwinding  in  all  three  direc¬ 
tions.  This  feature  coupled  with  a  simplified  left-hand-side  has  produced  an  efficient 
and  accurate  numerical  scheme  that  is  about  five  times  computationally  efficient.  These 
additional  features  are  utilized  from  the  works  of  Ref.  28.  The  use  of  Roe’s  upwinding 
eliminates  the  addition  of  explicit  numerical  dissipation  to  control  the  stability  of  the 
code  and  thus  is  known  to  produce  less  dissipative  solutions.  This  feature,  coupled  with 
a  fine  grid  description  in  the  tip  region  would  enable  better  description  of  the  wake  sim¬ 
ulation.  It  is  hoped  that  the  preservation  of  the  wake  structure  in  the  flowfield,  at  least 
in  the  initial  formation  and  roll-up  stages,  would  be  better  in  the  absence  of  additional 
numerical  viscosity. 

Body-conforming  finite-difference  grids  have  been  used  in  the  present  study. 
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The  helicopter  rotor  blades  considered  are  of  rectangular  plan  form  with  no  twist  or 
taper.  The  grid  topologies  consist  of  warped  spherical  0-0  for  the  first  problem  and 
C-H  cylindrical  grid  for  the  second  problem.  For  the  first  problem,  the  flowfield  grid 
is  generated  using  the  three-dimensional  hyperbolic  grid  solver  of  Steger  and  Chaussee 
(Ref.  29)  with  proper  clustering  in  the  leading  and  trailing  edge  regicrs  and  in  the 
tip  region.  The  grid  is  nearly  orthogonal  at  the  surface  and  the  spacing  in  the  normal 
direction  at  the  surface  is  chosen  to  be  0.0002C  for  the  viscous  grid  and  0.02C  for  the 
Euler  grid.  The  viscous  computations  are  done  on  a  grid  consisting  of  21  points  along  the 
span,  101  in  the  periodic  direction  around  the  blade  and  31  in  the  normal  direction.  The 
Euler  grid  is  more  coarser  with  only  21x101x15  grid  points.  The  outer  grid  boundary 
was  chosen  to  be  at  8  chords  in  all  directions  for  both  viscous  and  inviscid  grids. 

The  C-H  cylindrical  grid  used  for  the  lifting  hover  calculation  appears  to  be  the 
desirable  topology  because  of  the  cylindrical  nature  of  the  flow.  The  grid  is  clustered 
as  before  to  enable  better  resolution  of  the  flow.  Th  grid  used  has  approximately  one 
million  points  with  217  in  the  chordwise  direction,  71  in  the  spanwise  direction  and  61 
in  the  normal  direction.  The  spacing  near  the  wall  is  0.C0005.  The  smallest  spacing  at 
the  bevel  tip  in  the  spanwise  direction  is  of  the  order  of  the  semi-thickness  of  the  blade. 
To  test  the  code  a  smaller  grid  was  constructed  by  selecting  every  other  point  in  all  the 
three  directions. 

The  boundary  conditions,  both  surface  and  farfield,  are  applied  explicitly.  At  the 
rotating  blade  surface  because  of  the  noslip  condition  for  viscous  flows,  the  contravariant 
velocities  U ,  V  and  W  are  set  to  zero.  However,  the  surface  velocity  in  the  inertial  frame 
is  nonzero  since  it  is  equal  to  the  surface  grid  velocity,  determined  by  the  time  metrics 
£t,  r/i  and  as  the  blade  (and  the  grid  attached  to  it)  moves  in  azimuth.  The  Euler 
calculations  use  an  extrapolation  of  contravariant  velocities  to  the  surface.  The  density 
at  the  wall  is  determined  by  zeroth-order  extrapolation.  The  pressure  along  the  body 


surface  is  calculated  from  the  normal  momentum  relation  (see,  for  example  Ref.  22). 
Having  calculated  the  density  and  pressure,  the  total  energy  is  determined  from  the 
equation  of  state. 

At  the  farfield  boundary  the  flow  quantities  are  either  fixed  or  extrapolated  from 
the  interior  depending  on  whether  the  flow  is  subsonic  or  supersonic  and  if  it  is  of  inflow- 
or  outflow-type  at  the  boundary.  The  characteristic  velocities  of  the  Euler  equations 
determine  the  number  of  flow  properties  to  be  specified  to  control  the  reflections  of  waves 
from  the  boundaries.  For  the  subsonic-inflow  boundary,  four  quantities  must  be  specified 
and  one  quantity  is  determined.  The  four  specified  here  are  a  Riemann  invariant,  the 
entropy,  and  two  tangential  velocities;  the  quantity  that  is  calculated  is  also  a  Riemann 
invariant.  For  the  supersonic-inflow,  all  flow  quantities  are  specified.  At  the  subsonic- 
outflow  boundaries,  only  one  quantity  is  specified.  For  the  supersonic-outflow  condition 
all  flow  quantities  are  extrapolated  from  the  interior.  The  plane  containing  the  blade 
root  is  chosen  very  close  to  the  rotation  axis  of  the  blade  (at  R  =  1.0(7)  and  is  also 
treated  as  a  farfield  boundary. 

For  the  hovering  rotor  problem,  the  farfield  boundary  is  updated  explicitly  by 
means  of  characteristic-type  boundary  conditions  procedure  (Ref.  30),  although  the 
effect  of  using  this  procedure  was  found  to  be  unimportant.  This  is  because  the  Roe’s 
upwinding  procedure  used  in  the  numerical  scheme  would  otherwise  treat  the  boundaries 
in  a  1-D  characteristic  sense.  At  the  inboard  plane,  near  the  center  of  rotation,  a  linear 
extrapolation  is  used  for  density  and  the  three-velocities,  and  the  energy  is  calculated 
from  the  equation  of  state  by  prescribing  the  pressure.  A  periodicity  condition  is  used 
to  swap  the  information  at  the  front  and  back  boundaries  of  the  cylidrical  grid  topology 
(  for  a  two  bladed  rotor)  for  wake  capture  feature.  At  the  bottom  boundary  of  the 
cylindrical  grid,  a  boundary  condition  similar  to  that  prescribed  at  the  farfield  boundary 
in  the  first  problem  is  imposed  explicitly. 


11 


The  interacting  line  vortex  generated  upstream  of  the  rotor  by  a  rectangular 
wing  is  fixed  in  the  inertial  space  along  a  line  of  constant  x.  It  is  assumed  to  have  an 
analytical  representation  for  the  cylindrical  velocity  distribution  given  by 


ve  = 


CkT 

2irr 


(to) 


(10) 


The  constant  Ck  in  the  above  expression  has  been  determined  to  be  equal  to  0.8  by 
matching  the  peak  tangential  velocity  with  the  experimentally  measured  value  at  the 
measured  radial  distance  reported  in  Ref.  6.  The  constant  a0  is  approximately  equal 
to  the  radius  of  the  viscous  core  of  the  vortex  (equal  to  0.167C  in  the  present  case). 
The  induced  velocity  due  to  this  line  vortex  is  calculated  using  the  Biot-Savart  Law  and 
the  pressure  field  is  calculated  by  solving  the  radial  momentum  equation;  the  vortex  is 
initialized  well  upstream  of  the  interaction  region  as  described  in  Ref.  5. 
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3.  RESULTS  AND  DISCUSSION 


All  the  calculations  performed  in  the  first  part  of  this  study  are  done  in  a  time- 
accurate  fashion.  The  test  cases  calculated  correspond  to  the  two  sets  of  experiments  of 
Caradonna  et  al.  (Refs.  18-19)  for  parallel  and  oblique  blade-vortex  interactions.  Since 
one  of  the  purposes  of  the  experiments  was  to  collect  data  to  validate  numerical  methods, 
the  experimental  apparatus  was  kept  simple  to  ease  the  representation  by  numerical 
methods.  The  rotor  geometry  consisted  of  a  two-bladed  rigid  rotor  of  approximately 
14-chord  diameter.  The  blades  have  a  rectangular  planform  and  consist  of  an  NACA 
0012  airfoil  with  no  twist  or  taper.  The  rotor  blades  are  essentially  nonlifting  in  the 
absence  of  the  vortex  interaction.  The  interacting  vortex  was  generated  upstream  of  this 
rotor  by  a  lifting  rectangular  wing  with  an  NACA  0015  airfoil.  Parallel  and  oblique  BVI 
are  generated  by  positioning  the  line  vortex  along  the  y-axis  (z0  =  0.0)  and  ahead  of  the 
y-axis  (zQ  <  0.0)  so  that  the  rotor  blade  encounters  the  vortex  on  the  advancing  side  in 
the  first-  and  second-quadrants  of  the  azimuthal  travel  as  shown  in  Fig.  2. 

As  mentioned  before,  the  Navier- Stokes  hover  code  described  in  Ref.  21  provides 
the  starting  point  for  the  present  investigation.  In  this  study,  this  code  has  been  validated 
against  the  experimental  data  of  Caradonna  and  Tung  (Ref.  20)  for  the  nonlifting  and 
lifting  hover  conditions  for  both  the  subcritical  and  the  supercritical  tip  speeds.  Although 
the  main  focus  of  the  present  investigation  is  to  calculate  both  parallel  and  oblique  blade- 
vortex  interactions,  it  is  necessary  to  first  calculate  the  base-line  rotor-alone  flowfield 
solutions  (in  the  absence  of  vortex  interaction)  at  the  same  free  stream  conditions.  Two 
sets  of  results  for  parallel  blade- vortex  interactions,  corresponding  to  Mup  =  0.6  and  0.8 
and  /x  =  0.2  and  one  set  of  results  for  oblique  blade- vortex  interaction  corresponding  to 
Aftip  =  0.763  and  /x  =  0.197,  will  be  presented  in  the  following  sections  and  compared 
with  experimental  data.  At  the  Reynolds  numbers  corresponding  to  the  blade  tip  speeds 
in  these  experiments,  the  boundary  layer  on  the  entire  rotor  blade  could  be  assumed  to 
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be  turbulent  and,  as  mentioned  before,  an  algebraic  eddy  viscosity  model  of  Baldwin 
and  Lomax  (Ref.  23)  is  used  for  calculating  the  turbulent  viscosity. 

3.1  Rotor-Alone  Case 

The  calculation  of  the  rotor-alone  flowfield  solutions  serves  two  purposes.  First, 
it  enables  an  understanding  of  the  importance  of  unsteady  time-lag  effects  in  shock  wave 
growth  and  decay  as  well  as  the  three-dimensionality  of  the  flowfield  of  the  advancing 
rotor,  particularly  when  the  flow  is  supercritical.  Second,  it  provides  the  baseline  solution 
to  eventually  calculate  the  vortex  interaction  flowfield. 

Figure  3  shows  the  instantaneous  surface  pressure  distributions  for  various  az¬ 
imuthal  positions  of  the  rotor  for  the  case  of  Map  =  0.6  and  /x  =  0.2.  Examination  of 
the  results  suggests  that  at  this  reference  station  of  tb  =  0.893  and  for  the  subcritical 
flow  condition,  the  flow  behaves  as  if  it  is  quasi-steady  and  quasi-two-dimensional.  In 
spite  of  the  gradual  increase  of  the  blade-element  relative  speed  in  the  first  quadrant  and 
gradual  decrease  of  the  same  in  the  second  quadrant,  the  flowfield  appears  to  remain 
nearly  the  same  at  all  azimuth  locations  (as  seen  in  Fig.  3b)  indicating  that  the  un¬ 
steady  time-lag  effects  are  virtually  absent  for  this  flow  condition  and  the  flow  behaves 
as  if  it  is  quasi-steady.  In  fact,  also  shown  in  Fig.  3b  is  the  quasi-steady  solution  for  one 
azimuth  location  of  V’  =  90.°  which  is  in  perfect  agreement  with  the  unsteady  solutions 
at  different  azimuthal  positions.  Also,  the  comparison  of  three-dimensional  results  of 
Fig.  3b  with  the  two-dimensional  results  (Ref.  6)  of  Fig.  3a  shows  that  these  results  are 
nearly  identical  suggesting  that  the  flow  also  behaves  as  if  it  is  quasi-two-dimensional. 

At  the  supercritical  tip  flow  conditions,  corresponding  to  Mtip  =  0.8  and  /i  = 
0.2,  the  basic  rotor  flowfield  is  dominated  by  the  presence  of  a  strong  shock  wave  on  the 
advancing  blade  over  large  parts  of  the  first-  and  second-quadrants.  For  this  advancing 
rotor,  the  Euler  and  Navier-Stokes  results  of  instantaneous  surface  pressures  are  shown 
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in  Fig.  4  at  the  radial  station  rs  —  0.893  for  different  azimuthal  positions  of  the  blade. 
The  Euler  results  show  good  agreement  with  the  experiments  and  the  Navier-Stokes 
results  have  only  a  qualitative  agreement  because  of  the  coarse  grid.  It  appears  that 
the  shock  wave  position  is  under  predicted  for  the  Navier-Stokes  method  and  slightly 
over  predicted  for  the  Euler  method.  For  this  coaxse  grid  solution,  the  Navier-Stokes 
method  predicts  the  decay  of  the  shock  much  faster  in  the  second-quadrant.  On  the 
other  hand,  the  Euler  method  predicts  that  the  shock  wave  persists  at  the  ip  =  180° 
azimuth.  Previous  two-dimensional  Navier-Stokes  calculations  of  the  same  flow  (Ref. 
6)  had  shown  strong  three-dimensionality  and  unsteady  time-lags  in  shock  wave  growth 
and  decay.  The  two-dimensional  assumption  for  this  flow  essentially  over-predicted  the 
shock  wave  position  and  strength,  unlike  the  subcritical  flow  condition.  In  contrast,  the 
three-dimensional  Euler  results  seem  to  follow  the  experimental  observation  correctly. 

Thus  at  this  flow  condition,  the  flowfield  is  highly  three- dimensioned  exhibiting 
strong  unsteady  time-lag  behaviour  in  the  shock  wave  formation  and  eventual  demise. 
The  unsteady  time-lag  feature  of  the  flowfield  is  demonstrated  in  the  surface  pressure 
distributions  of  Fig.  4  shown  for  ip  =  60°  and  ip  =  120°  azimuthal  positions,  which 
shows  that  the  two  results  are  very  different  from  each  other. 

In  principle,  the  shock  wave  should  attain  its  maximum  strength  when  the  rela¬ 
tive  flow  speed  reaches  a  maximum  value  for  ip  =  90°  blade  azimuthal  position.  However, 
the  shock  wave  keeps  growing  in  strength  even  after  the  blade  has  passed  the  ip  —  90°  az¬ 
imuthal  position  and  a  weak  shock  wave  seems  to  persist  even  at  the  ip  =  180°  azimuth. 
The  strengthening  of  the  shock  wave  in  the  first-quadrant  and  beyond,  and  the  slow 
demise  of  it  in  the  second-quadrant,  suggests  the  existence  of  a  strong  unsteady  time-lag 
effect  at  this  flow  condition.  Such  a  behavior  was  totally  absent  at  the  subcritical  flow 
condition.  The  presence  of  the  shock  wave  seems  to  introduce  not  only  the  time-lag 
in  the  adjustment  of  the  flow  as  the  blade  sweeps  in  azimuth,  but  also  makes  the  flow 
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highly  three-dimensional. 


In  contrast  to  the  flow  on  the  advancing  side  of  the  rotor  blade,  the  flowfield  of 
the  retreating  side  appears  really  benign  for  this  nonlifting  rotor.  The  lingering  effects  of 
the  shock  wave,  persisting  at  xp  =  180,°  soon  die  out  as  the  blade  sweeps  into  the  third- 
quadrant.  Shown  in  Fig.  5  are  the  surface  pressure  distributions  for  several  azimuthal 
locations  of  the  blade  in  third-  and  fourth-quadrants  on  the  retreating  cycle.  Since  the 
flowfield  is  basically  subcritical  on  the  retreating  side,  the  plots  of  Fig.  5  appear  to  have 
the  same  behavior  as  those  of  Fig.  3. 

Figure  6  shows  the  surface  pressure  distributions  on  the  advancing  rotor  at  rg 
=  0.946  for  the  flow  conditions  of  Mtiv  =  0.763  and  /x  =  0.197.  Since  this  flow  condition 
is  similar  to  that  of  the  case  of  Fig.  4,  it  has  the  same  qualitative  behaviour  of  the 
flow  on  the  advancing  side.  As  seen  in  these  plots,  the  Euler  calculations  are  in  good 
agreement  with  the  experimental  data.3  As  before,  the  Navier-Stokes  results  have  only 
a  qualitative  agreement  for  the  reasons  mentioned  before.  It  is  important  to  mention 
that  the  rotor-alone  case  involves  a  nonlifting  symmetrical  blade,  for  which  the  pressure 
distributions  should  be  identical  on  the  upper  and  lower  surfaces.  However,  the  data  in 
Fig.  6,  from  Ref.  19,  seem  to  have  more  scatter  than  the  corresponding  experimental 
results  from  Ref.  18. 

3.2  Quasi-steady  vs  Unsteady  Results 

Although  the  flowfield  of  a  lifting  helicopter  rotor  is,  in  general,  unsteady  and 
three-dimensional,  a  close  inspection  of  the  subsonic  surface  pressure  results  presented  in 
Fig.  3b  suggests  that  the  unsteady  and  quasi-steady  results  are  nearly  the  same  for  this 
nonlifting  flow  condition.  In  fact,  a  quasi-steady  flow  calculated  by  freezing  the  rotor 
blade  at  certain  azimuthal  locations  (say,  xp  =  0.0°  )  produces  nearly  the  same  flowfield 
as  an  unsteady  flowfield  irrespective  of  the  flow  conditions.  In  the  present  approach, 
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the  unsteady  calculations  are  started  from  a  baseline  quasi-steady  solution.  At  the 
subcritical  flow  conditions,  similar  to  that  presented  in  Fig.  3,  an  unsteady  calculation 
can  be  started  from  any  azimuthal  location  of  the  blade  because  of  the  quasi-steady  and 
the  unsteady  solutions  being  nearly  identical  at  all  azimuthal  locations.  However,  this  is 
not  true  for  supercritical  flow  conditions  where  shock  waves  are  present  in  the  flowfield. 
For  example,  Fig.  7  shows  the  quasi-steady  and  unsteady  surface  pressure  distributions 
for  the  case  of  Mtip  =  0.8  and  [l  =  0.2  at  azimuthal  locations  of  ip  =  60°  and  120°. 
For  this  flow  condition,  these  results  indicate  that  the  quasi-steady  assumption  over 
predicts  the  shock  strength  at  rj>  =  60°  and  under  predicts  at  ip  =  120°.  Nevertheless, 
at  =  0.0°  position,  the  quasi-steady  and  unsteady  results  are  nearly  identical  even  for 
this  flow  condition.  So,  for  flow  conditions  that  are  severe,  it  is  still  reasonable  to  start 
the  unsteady  marching  from  the  quasi-steady  solution  calculated  at  ip  =  0.0°  azimuthal 
location. 

3.3  Parallel  Blade- Vortex  Interaction 

During  the  unsteady  three-dimensional  close  encounter  of  a  curved  tip  vortex 
with  a  rotating  blade,  the  helicopter  rotor  undergoes  a  variety  of  blade-vortex  inter¬ 
actions  depending  on  the  interaction  angle  between  the  leading  edge  of  the  rotating 
blade  and  the  curved  line  vortex.  These  interactions  are  generally  unsteady  and  three- 
dimensional.  One  limiting  case  of  such  an  encounter,  for  zero  interaction  angle,  is  termed 
parallel  interaction  (see,  for  example,  Ref.  5).  In  the  experimental  configuration  consid¬ 
ered  here,  this  interaction  occurs  around  rj>  =  180°  azimuthal  position.  For  an  observer 
riding  with  the  blade  (at  a  given  reference  station  along  the  span),  it  appears  as  though 
the  observer  is  passing  a  fixed  vortex  in  the  flow.  For  this  reason,  this  interaction  is 
sometimes  approximated  as  two-dimensional  and  unsteady. 

To  calculate  accurately  the  blade-vortex  interaction  flowfield,  it  is  necessary  to 
preserve  the  vortex  structure  without  numerical  diffusion.  As  mentioned  before,  one 
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method  that  has  been  demonstrated  to  work  effectively  and  economically  in  achieving 
this  is  the  prescribed-disturbance  scheme  (Ref.  5).  The  effectiveness  of  this  scheme  is 
illustrated  in  Fig.  8  (reproduced  from  Ref.  7)  showing  the  variation  of  lift  coefficient  as 
a  function  of  the  vortex  location  during  an  airfoil- vortex  interaction.  Also  shown  in  this 
figure  is  the  Euler  lift  distribution  calculated  using  a  conventional  (vortex  capturing) 
technique  (Ref.  31).  Although  the  two  Euler  solutions  are  computed  on  the  same  grid, 
the  numerical  dissipation  associated  with  the  finite  grid-spacing  progressively  weakens 
the  gradients  and  reduces  the  effective  vortex  strength  in  the  conventional  method.  Also, 
this  numerical  error  is  grid-dependent  (the  finer  the  grid  the  lesser  the  error);  however,  it 
is  completely  absent  in  the  prescribed- vortex  solutions,  which  are  essentially  independent 
of  the  grid.  Hence,  the  prescribed-disturbance  method  has  been  used  here  for  preserving 
the  vortex  structure. 

3.3.1  Subcritical  Interaction 

The  results  of  a  subcritical  parallel  BVI  are  discussed  here.  This  case  corre- 
sponds  to  the  flow  conditions  of  Mup  =  0.6,  n  =  0.2,  and  CfcT  =  0.133  at  a  blade 
reference  station  (marked  AA  in  Fig. 2)  tb  —  0.893.  The  interacting  vortex  is  located  at 
z0  =  0.0  and  z0  =  -0.4  along  the  y-axis.  To  calculate  the  BVI  flowfield,  the  interacting 
vortex  is  initialized  with  the  baseline  rotor  solution  at  the  azimuthal  position  of  ip  = 
120.°  The  unsteady  flow  field  is  monitored  as  the  blade  continues  to  advance  in  azimuth, 
and  the  peak  interaction  effects  are  shown  in  Fig.  9  in  the  form  of  instantaneous  surface 
pressure  distributions  measured  around  ip  =  180°  azimuthal  position. 

The  present  Euler  results  are  in  very  good  agreement  with  both  the  experimental 
data  (Ref.  18)  and  a  previous  two-dimensional  Navier-Stokes  calculation  (Ref.  6).  As 
seen  from  Fig.  9a  for  ip  =  178.15,°  the  lift  on  the  blade,  which  is  initially  zero  in  the 
absence  of  vortex  interaction,  is  negative  (because  of  the  sense  of  rotation  of  the  vortex) 
and  gradually  becomes  positive  as  the  blade  passes  the  line  vortex.  This  crossover  of  lift, 
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from  negative  to  positive,  seems  to  occur  when  the  vortex  is  approximately  aligned  with 
the  quarter-chord  line  of  the  blade  (corresponding  to  xv  =  0.0).  As  mentioned  before, 
the  two-dimensional  approximation  of  this  three-dimensional,  unsteady  interaction  is,  in 
fact,  a  very  good  assumption  for  this  subcritical  flow  condition.  The  close  agreement  of 
the  viscous  and  inviscid  results  suggests  that  the  viscous  effects  are  unimportant  for  this 
relatively  weak  interaction. 

3.3.2  Supercritical  Interaction 

>*s. 

This  case  corresponds  to  the  flow  condition  of  Mtip  =  0.8,  /z  =  0.2,  and  CfcT 
=  0.177  at  a  blade  reference  station  of  tb  —  0.893.  The  interacting  vortex  location 
is  same  as  in  subcritical  case.  This  flow  condition,  in  constrast  to  the  subcritical  flow 
condition,  exhibits  strong  unsteady  time-lag  effects  in  shock  wave  growth  and  decay  for 
the  advancing  rotor.  This  feature  seems  to  accentuate  the  three-dimensional  nature  of 
the  flow  and  a  two-dimensional  approximation  of  this  flow  over  predicts  the  shock  wave 
strength  and  location  as  demonstrated  in  Ref.  6. 

As  before,  to  calculate  the  BVI  flowfield,  the  interacting  vortex  is  initialized 
at  about  ip  =  120°  azimuthal  position  and  the  time  evolution  of  the  flow  is  monitored 
as  the  blade  advances  in  azimuth.  Figures  10-11  show  the  interaction  flowfield  results 
in  terms  of  instantaneous  surface  pressures  for  several  azimuthal  positions  of  the  blade 
for  the  Euler  and  the  Navier-Stokes  calculations.  Also  shown  in  these  figures  are  the 
experimental  data  (Ref.  18).  Comparison  of  the  results  shows  that  while  the  Euler  results 
(Fig.  10)  have  good  agreement  with  the  experimental  data  for  all  azimuthal  locations  of 
the  blade,  the  Navier-Stokes  results  (Fig.  11)  show  only  a  qualitative  agreement  with  the 
experimental  data.  This  feature  is  strongly  dependent  on  the  accuracy  of  the  base-line 
solution  at  the  azimuthal  location  where  the  vortex  is  initialized.  The  vortex  is  intialized 
at  the  ip  =  120°  position  of  the  blade  for  this  case.  Comparison  of  surface  pressures  of 
Fig.  4d  at  this  location  for  the  Navier-Stokes  and  Euler  calculations  shows  that  the 
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Euler  results  have  better  agreement  with  experiments  than  the  Navier-Stokes  results. 
Although  the  Euler  grid  is  very  coarse,  it  seems  to  capture  the  essential  details  of  the 
flow  sharply  including  the  shock  waves.  Since  the  Navier-Stokes  grid  is  too  coarse,  the 
shocks  appear  smeared.  Nevertheless,  it  captures  the  essential  qualitative  features  of  the 
flow. 


The  effect  of  vortex  interaction  is  to  induce  time  dependent  aerodynamic  forces 
on  the  rotor.  For  example,  as  seen  from  the  surface  pressure  plots  of  Figs.  10-11,  the  lift 
on  the  blade,  which  is  intially  zero,  becomes  negative  (because  of  the  sense  of  rotation 
of  vortex)  due  to  induced  downwash  and  gradually  increases  to  a  positive  value  as  the 
blade  passes  the  fixed  vortex.  The  peak  effects  of  the  interaction  appear  to  occur  when 
the  blade  leading  edge  is  approximately  above  the  line  vortex. 

3.4  Oblique  Blade- Vortex  Interaction 

This  case  cc~responds  to  a  free  stream  condition  of  Mtip  =  0.763  and  \x  = 
0.197.  The  interacting  vortex  located  at  x0  =  -2.13  and  z0  =  -0.25  has  a  strength  of 
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CkT  =  0.179  for  the  reference  blade  section  at  =  0.946.  In  contrast  to  the  parallel 
interaction,  the  oblique  interaction  occurs  over  a  larger  azimuthal  sweep  of  the  rotor 
blade  starting  in  first-quadrant  and  completing  in  second-quadrant.  The  peak  interaction 
occurs  at  around  ip  —  20°  in  the  first-quadrant  and  around  rp  =  160°  in  the  second- 
quadrant.  Therefore,  the  vortex  is  initialized  at  tp  =  0.0°  and  the  evolution  of  the 
unsteady  interacting  flow  is  monitored  for  the  advancing  rotor. 

Figures  12-13  show  the  surface  pressure  distributions  for  several  azimuthal  po¬ 
sitions  of  the  blade  from  the  Euler  and  the  Navier-Stokes  calculations  of  the  interaction 
flowfield.  As  observed  before,  the  transonic  Euler  calculations  agree  better  with  the 
experimental  data3  than  the  Navier-Stokes  calculations  for  the  reasons  cited  earlier.  On 
the  other  hand,  the  interaction  in  the  first-quadrant  occurs  at  subcritical  flow  condition 
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and  therefore  the  agreement  of  both  the  Euler  and  Navier-Stokes  calculations  with  the 
experimental  data  is  good  at  ip  =  50°.  For  the  interaction  in  the  sccond-qudrant,  where 
the  flow  is  supercritical,  the  Euler  results  have  better  agreement  with  the  experiments 
than  the  Navier-Stokes  results.  Compared  to  the  preceding  transonic  parallel  interac¬ 
tion  results,  the  Navier-Stokes  results  for  this  interaction  appear  to  agree  better  with 
the  experimental  data.  This  is  in  direct  consequence  of  the  better  agreement  of  the 
rotor-alone  Euler  and  the  Navier-Stokes  surface  pressures  at  all  azimuthal  positions,  as 
seen  in  Fig.  12  for  this  flow  condition.  The  maximum  relative  blade  tip  speed  in  this 
case  (for  ip  =  90°  position)  is  0.91  compared  to  a  value  of  0.96  corresponding  to  the 
transonic  parallel  interaction  case  and  is  thus  a  weaker  interaction. 

From  the  surface  pressures  of  Fig.  13,  it  can  be  seen  that  the  lift  on  rotor  blade, 
which  is  initially  zero,  becomes  negative  due  to  the  vortex  influence  (because  of  the  sense 
of  rotation),  increases  to  a  maximum  (negative)  value  in  second-quadrant,  and  changes 
over  to  a  positive  value  at  about  ip  =  165°  azimuthal  position  of  the  blade.  The  peak 
interaction  appears  to  occur  around  the  ip  =  160°  azimuthal  position. 

3.5  Lifting  Hovering  Rotor  Flowfield 

This  section  presents  the  steady  solutions  of  a  lifting  hovering  rotor.  The  flow- 
field  of  a  hovering  rotor  can  be  either  steady  or  unsteady  depending  on  the  reference 
frame  of  the  observer.  It  is  steady  if  the  observer  is  fixed  to  the  blade  and  riding  with 
it.  On  the  other  hand,  if  the  observer  is  stationed  in  the  inertial  reference  frame,  the 
flowfield  appears  unsteady.  To  take  advantage  of  the  steady  state  of  the  flowfield,  a  vari¬ 
able  time-step  option  (see  Ref.  31)  is  used  to  accelerate  the  convergence  criteria  of  the 
numerical  system.  Rest  of  the  features  of  the  physical  problem  including  the  unsteady 
metrics  are  identical  to  the  first  problem  presented  above.  The  major  advantage  in  this 
steady  state  problem  is  that  the  blade  need  not  be  physically  advanced  at  each  time 
step 
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The  test  cases  selected  for  calculations  correspond  to  the  experimental  test  con¬ 
ditions  of  a  two-bladed  rigid  rotor  (Ref.  20).  The  rectangular  blades  of  this  rotor  are 
made  up  of  NACA  0012  airfoil  sections  and  have  an  aspect  ratio  of  6.  The  blades  are 
untwisted  and  untapered  as  before.  Three  typical  sets  of  conditions  are  considered  rep¬ 
resenting  subcritical  and  supercritical  tip  speeds  with  two  different  collective  pitches. 
The  Reynolds  number  variation  in  these  sets  range  from  1.92  to  3.83  million  based  on 
the  tip  speeds.  In  the  calculation  procedure,  the  flowfield  is  initially  quiescent  and  the 
evolution  of  the  flowfield  is  monitored  as  the  blade  is  set  in  motion. 

Figures  14-18  show  steady  surface  pressure  distributions  for  the  subcritical  case 
of  Mtip  =  0.44  with  a  collective  pitch  of  8°  and  a  flow  Reynolds  number  based  on  the 
tip  speed  of  1.92x10®  at  various  radial  stations  along  the  blade.  The  solid  line  which 
is  the  present  calculations  shows  excellent  agreement  with  the  experimental  data  (Ref. 
20)  at  all  radial  stations.  Also  shown  in  these  is  a  previous  Navier-Stokes  calculation 
of  Srinivasan  and  McCroskey  (Ref.  21)  that  used  a  simple  wake  model  to  bring-in  the 
influence  on  the  wake.  Although  there  is  good  qualitative  agreement  between  these  two, 
the  present  calculation,  which  does  not  use  any  ad  hoc  assumptions  on  the  wake  and  its 
influence,  are  in  very  good  agreement  with  the  experimental  data.  The  two  calculations 
are  done  on  a  comparable  size  grids  but  have  different  grid  topologies. 

The  present  wake  capture  calculations  are  atleast  5  times  cheaper  in  computer 
time  requirements  compared  to  the  one  presented  in  Ref.  21.  Calculations  have  also  been 
done  on  a  smaller  grid,  which  was  1/8  the  size  of  the  present  grid  of  about  a  million  point, 
and  the  results  show  very  good  agreement  with  experiments  even  on  this.  The  deficiency 
is  in  the  inboard  part  of  the  blade  where  the  grid  is  not  very  fine.  The  calculations  on 
the  smaller  grid  takes  typically  less  than  one  hour  on  Cray2  supercomputer.  This  time 
is  comparable  to  even  many  potential  flow  codes  (Refs.  32-34)  that  are  in  use  today. 

Figures  19-22  show  steady  surface  pressures  for  the  supercritical  case  of  Mtip 
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=  0.877  with  a  collective  pitch  of  8°  and  corresponding  Reynolds  number  of  3.83x10®. 
The  surface  pressures  are  shown  for  four  radial  stations  and  are  compared  with  the 
experimental  data  of  Caradonna  and  Tung  (Ref.  20)  and  the  Navier-Stokes  calculations 
of  Srinivasan  and  McCroskey  (Ref.  21).  As  seen  in  these  figures,  the  present  wake  capture 
Navier-Stokes  results  are  in  much  better  agreement  with  the  experimental  data  than  the 
Navier-Stokes  results  of  Ref.  21.  These  results  further  demonstrate  the  superiority  of 
this  present  methodology  in  predicting  accurate  results  with  out  any  ad  hoc  inputs  to 
describe  the  wake  which  is  an  important  influence  in  the  lifting  hovering  rotor  flowfields. 

Figures  23-27  show  surface  pressures  for  the  conditions  of  Mt;p  =  0.877,  collective 
pitch  of  12°,  and  flow  Reynolds  of  3.55x10®.  These  calculations  are  done  on  the  coarser 
grid.  The  results  show  only  reasonable  agreement  with  the  experimental  data.  This  is 
a  severe  case  in  terms  of  vorticity  in  the  wake.  These  results  are  in  similar  agreement 
with  experimental  data  as  some  of  the  potential  flow  methods  that  use  wake  models 
(Ref.  33-34).  It  should  be  emphasised  that  the  present  calculations,  with  the  smaller 
grid,  takes  less  than  one  hour  of  computer  time  on  a  supercomputer  and  no  additional 
inputs  are  made  to  describe  the  wake. 
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4.  CONCLUDING  REMARKS 


The  research  presented  here  address  two  related  problems  of  a  helicopter  ro¬ 
tor  and  are  therefore  three-dimensional,  unsteady  and  viscous  dominated.  For  the  first 
problem  considered  here,  a  numerical  procedure  is  presented  to  calculate  the  unsteady, 
viscous  and  inviscid,  three-dimensional  flowfields  of  a  helicopter  rotor  in  forward  flight 
and  parallel  and  oblique  blade-vortex  interactions  in  subcritical  and  supercritical  flow 
conditions.  Important  flow  features  such  as  unsteady  time-lag  effects  in  shock  wave 
growth  and  evential  decay,  as  well  as  the  importance  of  three-dimensional  effects  are 
discussed.  While  it  is  possible  to  approximate  the  parallel  blade-vortex  interaction  as 
two-dimensional  and  unsteady  when  the  entire  flow  is  subcritical,  this  same  problem 
under  supercritical  flow  condition  and  also  the  oblique  blade- vortex  interaction  are  nat¬ 
urally  three-dimensional,  and  unsteady  and  has  to  be  treated  that  way  without  any 
simplifications. 

The  numerical  results  are  compared  with  two-sets  of  test  data  generated  by 
Caradonna  et  al.  (Refs.  18-19)  on  a  model  two-bladed  rigid  rotor  in  a  wind  tunnel.  While 
the  Euler  results  have  good  agreement  with  experiments,  the  Navier-Stokes  results,  with 
a  very  coarse  grid,  have  only  a  good  qualitative  agreement.  Nevertheless,  the  study 
demonstates  that  the  numerical  methodology  presented  here  is  capable  of  calculating 
accurately  the  flowfield  of  a  helicopter  rotor  including  the  effects  of  vortex  interaction 
under  subcritical  and  supercritical  flow  conditions. 

The  second  problem  considered  is  that  of  a  lifting  hovering  rotor  flowfield.  In 
contrast  to  the  currently  available  potential  flow  methods  that  use  wake  modeling  to 
do  these  flowfield  calculations,  the  Navier-Stokes  method  presented  here  uses  no  ad  hoc 
inputs  to  describe  the  wake.  The  wake  and  its  influence  are  captured  as  a  part  of  the 
overall  flow  field  solution.  The  method  is  very  fast  and  robust  and  is  applicable  to 
arbitrary  blade  planforms  and  geometries  with  out  any  modifications.  The  calculations 
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performed  in  this  report  for  subcritical  and  supercritical  flow  conditions  show  excellent 
agreement  with  the  experimental  data. 
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Fig.  1  Schematic  of  the  complex  flowfield  of  a  helicopter  rotor  in  forward  flight 
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Fig.  2  Schematic  of  the  rotor  BVI  wind  tunnel  test. 
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Fig.  3  Instantaneous  surface  pressure  distributions  at  the  reference  blade  section  for 
advancing  rotor.  Mtip  =  0.6,  fi  =  0.2,  Re  =  2.1x10®,  rB  =  0.893. 
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Fig.  5  Instantaneous  surface  pressure  distributions  at  the  reference  blade  section  for 
retreating-blade- Euler  solution.  Mup  =  0.8,  fj.  =  0.2,  tb  =  0.893. 
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-  NAVIER-STOKES  CALCULATIONS 


-  EULER  CALCULATIONS 

•  EXPERIMENT  (CARAOONNA  et  al.) 
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—  EULER  UNSTEADY 
EULER  QUASISTEADY 
•  EXPERIMENT  (CARAOONNA,  et  at.) 


Fig.  7  Comparison  of  pressure  distributions  for  quasi-steady  and  unsteady  Euler 
calculations.  Mtip  =  0.8,  n  =  0.2,  rB  =  0.893. 


—  EULER  (PERTURBATION! 
- EULER  (NONPERTURBATION) 


Fig.  8  Effectiveness  of  the  prescibed-vortex  method  for  two-dimensional  BVI: 
NACA  64A006  airfoil.  Af*  =  0.85,  f  =  0.2,  z0  =  -0.26. 
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- 3-0  EULER  RESULTS-PARALLEL  INTERACTION 

*  X  2-D  NAVIER  STOKES  RESULTS 

•  0  EXPERIMENT  (CARADONNA,  et  al.) 
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Fig.  9  Comparison  of  instantaneous  surface  pressures  during  subsonic  parallel  blade- 
vortex  interaction.  The  solid  line,  the  solid  circle,  and^the  symbol  +  refer  to  the  vortex 
side  of  the  blade.  Mtip  =  0.6,  n  =  0.2,  Re  =  2.1xl06,  f  =  0.133,  xe  =  0.0,  ze  =  -0.4, 
rg  =  0.893. 
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Fig.  10  Instantaneous  surface  pressures  during  transonic  parallel  blade-vortex  interaction 
-  the  Euler  results.  The  solid  line  and  the  solid  circle  refer  to  the  vortex  side  of  the  blade. 


Mtiv  =  0.8,  a  =  0.2,  f  =  0.177,  *o  =  0.0,  =  -0.4,  rB  =  0.893. 


X  X 

Fig.  11  Instantaneous  surface  pressures  during  parallel  blade- vortex  interaction  -  the 
Navier-Stokes  results.  The  solid  line  and  the  solid  circle  refer  to  the  vortex  side  of  the 
blade.  Mtip  =  0.8,  y  =  0.2,  Re  =  2.89xl06,  P  =  0.177,  x0  =  0.0,  zc  =  -0.4,  rB  =  0.893. 
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Fig.  12  Instantaneous  surface  pressures  during  oblique  blade- vortex  interaction  -  the 
Euler  results.  The  solid  line  and  the  solid  circle  refer  to  the  vortex  side  of  the  blade. 
Mtip  =  0.763,  n  =  0.197,  Re  =  2.75x10®,  f  =  0.179,  x„  =  -2.13,  z<>  =  -0.25,  rB  =  0.946. 
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Fig.  13  Instantaneous  surface  pressures  during  oblique  blade-vortex  interaction  -  the 
Navier-Stokes  results.  The  solid  line  and  the  solid  circle  refer  to  the  vortex  side  of  the 
blade.  Mtip  =  0.763,  y  =  0.197,  Re  =  2.75xl06,  f  =  0.179,  x0  =  -2.13,  z0  =  -0.25, 


tb  =  0.946. 
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Fig.  14  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Mtip  —  0.44, 
act  =  8°,  Re  =  1.92x10®  and  rg  =  0.5. 
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Fig.  15  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Mup  =  0.44, 
at  =  8°,  Re  =  1.92x10s  and  rB  =  0.68. 
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Fig.  16  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Mup  = 
at  =  8°,  Re  =  1.92x10*  and  r=  =  0.80. 
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Fig.  17  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Map  =  0.44, 
a«  =  8°,  Re  =  1.92x10®  and  rg  =  0.89. 
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Fig.  19  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Mtip  —  0.877, 
Oft  =  8°,  Re  =  3.83xl06  and  rB  =  0.5. 
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Fig.  20  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor. 
ott  =  8°,  Re  =  3.83x10®  and  rg  =  0.68. 
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Fig.  21  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Mtip  =  0.877, 
at  =  8°,  Re  =  3.83x10®  and  tb  =  0.80. 
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Fig.  22  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  M 
<Xt  =  8°,  Re  =  3.83x10s  and  rg  =  0.89. 
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Fig.  23  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Mu-p 
at  =  12°,  Re  =  3.53x10®  and  rg  =  0.5. 


$ 


9 


•  o  0.7940 
-  0.7940 


0.6800  Experiment-Caradonna  and  Tung 
0.6933  Navier-Stokes  Smallgrid  Results 


-1.2  -0.8  -0A  0.0 


•  o 


0.7940 

0.7940 


at  2y/b 

12.000  0.8900  Experiment-Caradonna  and  Tung 

12.000  0.8914  Navier-Stokes  Smallgrid  Results 


o 


Fig.  26  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Muv  =  0.794, 
at  =  12°,  Re  =  3.53x10®  and  rB  =  0.89. 
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Fig.  27  Comparison  of  surface  pressures  for  a  lifting  hovering  rotor.  Mi 
a‘  =  12°,  Re  =  3.53x10*  and  rn  = 
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